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abstract 

The focus of the study was to develop the micromechanical model associated with proper damage model to 
predict the overall mechanical behavior of fiber/matrix unidirectional laminates. The present and first 
investigation studies the influence of fiber -matrix interface on the behaviour of fiber reinforced composite lamina 
using micromechanical models. Mechanical properties El and E 2 are determined at various volume fractions. 

The second investigation studies the micro-thermo elastic behaviour of the square unit cell of a hybrid fiber 
reinforced composite lamina. Later this model is extended to predict the coefficients of thermal expansion of 
graphite -boron hybrid fiber reinforced lamina for various volume fractions. In the third investigation, an 
analytical solution of the thermal stresses for a fiber embedded in a matrix is presented based on the idea of the 
finite element and under some simplifying assumptions. The analytical solution to the problem is found for the 
case when the length of the embedded bar (fiber) is much greater than its radius, and the Young's modulus of the 
matrix is much less than that of the fiber. The problem is also solved numerically by means of finite element 
analysis using ANSYS 10.0. Both results are compared and it is shown that both approaches coincide very close 
qualitatively and quantitatively although significant discrepancies may appear at specific points for specific 
cases. For all above three cases 3-D finite element models have been developed from the representative volume 
elements of the composite which are in the form of square unit cells. The finite element software ANSYS 10.0 has 
been successfully executed to evaluate the properties. 
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I. INTRODUCTION 

There has been a considerable increase in the use of advanced composite materials in various industries 
in recent years. The reason for this increase can be attributed to the great improvement of the stiffness-to-weight 
ratio and strength-to- weight ratio in Composite materials. The development of composite materials with reduced 
weight and increased strength relative to the conventional metals or alloys, has played a critical role in achieving 
higher operating performance, long-life and reduced costs. One such composite system is fiber reinforced 
composite laminates. The properties of composites, e.g. strength and stiffness, are dependent on the volume 
fraction of the fibers, and the individual properties of the constituent materials the fiber and matrix. In addition, 
the variation of lay-up configurations of composite laminates allows the designer greater flexibility when 
incorporating composites into a structure. With this flexibility, however, comes complexity in analysis of 
composite structures. Particularly, the damage and failure progression in laminated composites is very 
complicated compared to that of conventional metallic materials. The stiffness and strength in composite 
structures may vary due to the damage initiation and propagation. The damage modes are matrix cracking, 
interface debonding and delamination during the loading. Consequently, the modeling of composite materials is 
more complex than that of traditional engineering materials. 

II. LITERATURE REVIEW 

The failure and the analysis of composite structures can be various. The failure modes can be investigated 
separately, but these modes generally arise together (e.g.: matrix cracking and delamination). Due to this 
experimental observation the main problem is to construct models that can follow more than one failure modes. 
To solve this problem the most proper method is the finite element method (FEM), because of its simplicity. The 
aim of this paper is to summarize a few publications about the damage and failure of composite materials and 
structures. The publications were grouped according to the mode of failure, the aim is to give a short description 
about the finite element modelling of the given failure mode. The models were reconstructed using the finite 
element code COSMOS/MWatching the structure of composites we can notice its heterogeneity. The composite 
material consists of fibers and matrix. The fibers are embedded into the matrix. During the manufacturing 
process the matrix is warmed up to high temperature, and then the fibers are added. Due to the high temperature 
and the mismatch between the thermal properties of the components there is interaction between the fibers 
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and the matrix. The formation of interface is a consequence of this operation. Hence the composite is a three- 
phased material. The length and the arrangement of the fibers are various. The most frequent ones are the long 
unidirectional continuous and the short unidirectional arrangements. Another important feature of the 
composites is the periodicity. This feature can make the FE analysis of the 

Mesoscale analysis: 

A scale between the macro- and microscale, which can be a micromodel that is not Periodic in all directions or, 
considering laminated structures, this may mean the analysis of one elementary ply. In fact there is another 
scale, named multiscale analysis, coupling the three ones above. In this case we need parallel computers with 
big capacity. Multiscale approach was presented by FEYEL and CHABOCHE in [4] for 
composites. 

2.1 Material Models for The Finite Element Analysis 

There are significant differences between the material properties of the fiber and matrix. The stiffness of fiber is 
much higher than matrix stiffness. Major part of the loading is transferred onto the fiber. Experimental datas 
prove that the global composite structure behaves as a linear elastic material [10]. If local analysis is applied, 
then the considering of heterogeneity is necessary. 


2.1.1 Materials Of The Reinforcing Fibers 

The material of fibers has a wide range of type 

(e.g.: glass, carbon... etc.). It is assumed that the glass fiber is isotropic, the carbon fiber is transversely isotropic 
linear elastic material. The mechanical and other properties of the fibers are available in the literature [10]. 


2.1.2 Material Modelling Of The Matrix 

The modelling of the matrix is a far complex problem. In general the matrix is nonlinear elastic isotropic and 
hysteretic material. In case of polymer matrix (e.g.: epoxy, vinyl ester, polyester... etc.) applying linear elastic 
material model can give wrong results during the FE analysis. Due to this experience it is useful to apply 
rheological material models for the matrix. These are the viscoelastic models that consist of the well-known 
Maxwell and Kelvin-Voight elements and the combinations of them [2,6]. Metal matrix composites (MMC) 
require using elastoplastic, viscoplastic and elasto viscoplastic material model [3,4]. Inelastic material models are 
reachable in [14]. Inelastic behaviour can be simulated by stiffness or strength degradation. Previous 
investigations describe that using linear material model damage occurs in higher loading than in the case of 
nonlinear model [19]. Because of that the nonlinear material law describes more accurately the stress -strain 
relationship than the assumed linear. Commercial FE programs usually contain nonlinear material models, but 
there is a chance to build in another [19]. For that a user defined subroutine and the modification of it is needed 
[15]. 


2.2.2 Microscale finite element models 

Unit cells and RVE (Representative Volume Element) models are the most convenient forms in the composite 
materials simple. In literature there are three basic types of analysis [11, 19]. 

Microscale Analysis: 

Modelling of the local constituents (fiber, matrix, interface) using the periodicity and the symmetry. 

Macroscale analysis: 

Modelling of the global composite structure using the homogenized material properties of the composite. 
Micromechanical analysis. Making unit cells depends on the geometry and the applied load. Depending on the 
structure and fiber orientation these models can be two- and three dimensional ones. Due to the heterogeneity 
different material models must be applied to describe the three-phased composite. The construction and the 
mesh size of the finite element model have a great effect on the results of the analysis. Constructing the FE 
model the prescription of proper boundary conditions is required. 



Fig. 2.1 2D finite element model, squm e filtoi' 
airtmgeinent 
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Generally the load of micromechanical models are displacement or strain controlled. This means that one of the 
boundary curve or surface of the model is loaded with displacement in a given direction [1,19]. Fig. 2.1 shows the 
making of 2D model, the cut is perpendicular to the longitudinal direction of the fibers. Assuming the periodicity 
of composite a double symmetrical unit cell can be made. Using symmetrical conditions the FE model is only 
quarter of the former. Fig. 2.1 also shows the proper boundary conditions. In Fig. 2.1 the fiber and the matrix are 
covered by triangular, the interface is covered by rectangular plane elements in plane stress state. Depending on 
the aim of the analysis, there are several cases when the cut must be longitudinal. 

If composite structure is built by ply’s with different orientation a 3D FE model is required [3, 19]. In Fig. 2. 2 a 
three-dimensional unit cell and finite element model are shown. Using the periodicity a 3D multi -cell can be 
made, the finite element model is just one eighth of the multi-cell. Proper boundary conditions can be seen also 
in Fig. 2. 2 Between fiber and matrix a perfect bonding assumed and no debonding occurs. Three of the surfaces 
cannot move at the normal direction. Remained surfaces can be loaded with given displacements. 



Fiber arrangement is various. Most of the FEM programs haven’t got a well-built in CAD system. These CAD 
systems do not guarantee smoothness and accuracy that is always required. Because of that a proper CAD system 
is applicable [16]. Most of the FE codes can handle models drawed in CAD. Complicated fiber arrangements 
(e.g.: woven fabric, basket weave and tri-axial braided unit cells) raise the claim to use CAD system. More 
information can be found in [16]. The former paper presents modelling of composite unit cells using CAD, and 
also describes construction of unit cells with logical (Boolean) operations. 


2.3 FAILURE ANALYSIS OF COMPOSITES 

Failure means that one of the components reaches the yield stress, then damage occurrence and progression can be 
observed. Damage can progress on different directions. Usually matrix cracking is the first damage process to take 
place since the matrix has the lowest stress to failure. A failure criterion is needed to establishing initial damage. This 
can be the maximum stress or strain criterion [19]. Damaged zones can propagate and reach the border of each other. 
Finite element tests have a wide application range. Most of the cases a simplified, loaded model is investigated by 
FEM, the distribution of stress and strain fields give enough information about the damage process. Models can be 
analyzed with holes and notches, considering stress concentration effect. The following sections give a short review 
about the finite element analysis of damage and failure processes in composites using previous publications. 


2.3.1 Modelling of interface damage 

Experimental observations show that the damage process can be initiated by the fiber-matrix interface. Damage 
of the interface causes further problems (e.g.: debonding, fiber pullout). The interface strength plays an important 
role in the damage process. Depending on the loading the interface failure can cause by the mismatch of the 
modulus and Poisson’s ratios. SCHULLER et al. presented a very useful method for modelling the interface 
debonding [12]. A single fiber is embedded in a necked matrix block, as shown in Fig. 2. 3 Mechanical loading 
can be uniaxial compression or tension. Compression test has two disadvantages. First the fiber is often breaks or 
buckles before interfacial debonding can occur. Another drawback is the additional transverse stress components, 
which makes the stress field non-axisymmetrical. Under compression the transverse component is tensile one 
between the fiber and the matrix. 
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Fig.2.3 Necked single-fiber composite specimens 
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Fig 2.4 Cylinder (1.), hexagonal (2.) and 
special composite (3.) models 
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Fig 2.4 Cylinder (1.), hexagonal (2.) and special composite (3.) models 


The hexagonal unit cell shows the results that neighboring fibers more or less influence the residual stress 
components. Results obtained by the composite model show a decrease of thermal residual stresses with 
increasing local fiber-volume fraction. The location and the value of applied stress at which failure occurs are 
influenced by the local fiber-volume fraction. Depending on the applied matrix criterion damaged area can be 
different. The three used criterions are the parabolic, the strain energy density and the von Mises criterion. 
Varying the fiber- volume fraction these criterions give very distinct results. The von Mises criterion predicts over 
the maximum load, so it is not suitable to predict initial matrix failure. In fact thermal residual stresses reduce the 
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maximum bearable load, but they do not affect the location of initial failure. 

The effect of thermal residual stresses in short-fiber-reinforced composites was investigated by ZHOU et al. [21] 
with Al matrix and SiC fiber. The applied FE model and the boundary conditions are shown in Fig.2.5. 

Fig 2.5 FE model for short-fiber-reinforced Al/SiC composite 

Cooling from high to room temperature the residual shear stress concentrates near the corner of the fiber end and 
induces initial damage along the interface. Considering thermal cycling the radial, longitudinal and shear stress 
distributions were investigated. The results show considerable stresses ahead of the fiber end and around the 
corner of the fiber. Applying thermal cycling interface damage is more dominant than matrix damage. Thermal 
cycling decreases the strength and the toughness of metal matrix composites (MMC). The paper shows that the 
thermal cycle range rather than the number of cycles dominate damage process. Cooling and thermal cycling 
change the mode of failure. Considering only tensile load failure arises as matrix damage around the fiber end 
and no debonding occurs. Considering cooling with tensile load debonding along the fiber end and transverse 
damage across the matrix occur. Finally applying thermal cycling and mechanical load the MMC fails by 
debonding along the fiber ends and fiber sides, followed by necking instability of the matrix ligament. 
Microscopic behaviour of MMC is strongly influenced by the inelastic deformation of the matrix. ISMAR and 
co workers presented a study about MMCs reinforced by different types of fiber [7]. The damage of matrix is 
affected by the properties of the fibers. A nonsymmetrical [0790°] fiber arrangement model was investigated, 
Fig. 2. 6 depicts the 3D model. Boundary conditions are identical with shown in Fig. 2. 2. The three different fiber 
types: a-A1203, SiC and carbon fiber. All fibers behave elastically, a-A1203 and SiC are isotropic, the carbon 
fiber is described as transversely isotropic material. Interface between fiber and matrix assumed to be very strong 
and no debonding occurs. The a-A1203 has the highest tensile modulus, SiC has the lower, and modulus of 
carbon fiber is between the formers. During fabrication process, the composite passes through several 
temperature regions. Thermal residual stresses cause local inelastic deformation. 
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Fig 2.6 3D FE model for metal matrix composite 


A cooling process is firstly examined. Seeing the three fiber types the inelastic strain shows highest values in the 
case of carbon fiber. After cooling the equivalent inelastic strain concentrates areas shown in Fig. 2. 6. Inelastic 
strain distribution reaches the maximum values where the fibers are the closest to each other. Varying the fiber- 
volume fraction does not significantly alter the relative distribution of equivalent inelastic strain but strongly 
affects the extent of it. Because of that the effects of thermal cooling considered the following investigations of 
the mechanical behaviour. Mechanical loading is displacement controlled with given strain rate. Under 
monotonic loading A1203 -reinforced composite shows maximum inelastic strain owing to the larger difference 
between the tensile modulus of the fiber and matrix. The stress-strain curve of the three fiber reinforced 
composites shows that the slowest decline of the tangent modulus can be observed in the case of carbon fiber. On 
the other hand A1203 fiber has the fastest decline. Cyclic loading was considered with given strain rate and 
temperature. Metal matrix behaves in elastically applying given loading amplitude. The density of inelastic 
performance after forty cycle decreases depending on the fiber type. Increasing loading amplitude inelastic 
behaviour strongly increases. Cyclic loading indicates damage progression and damage accumulation. Damage is 
distinctly concentrated in two matrix areas; the A1203 -reinforced composite exhibits the highest damage values. 
Heightening the fiber- volume fraction, the inelastic deformation slightly increases. Fiber with large tensile 
modulus increases the strength of composite. On the other hand the danger of premature failure during cyclic 
loading has to be considered. Finally the fiber-volume fraction decisively affects the mechanical behaviour of 
metal matrix composites. 
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2.3.3 Matrix damage and matrix cracking 

Matrix crack occurs when the equivalent stress reaches the yield stress, after fiber end can debonded and a local 
plastic zone in the matrix is the result of further loading. In the case of short-fiber reinforced composites around 
the fiber end a penny-shaped crack can raise. SIRIVEDIN et al. presented a study about matrix cracking in short- 
fiber composite [13]. Depending on the strength of the interface this crack can progress on various ways. In the 
case of relatively weak interface, a cylindrical interface crack propagates from the debonded fiber end. If the 
interface is strong then a conical matrix crack can be observed propagating from the debonded fiber end. This 
can be described by the angle of the conical crack. The applied finite element model is depicted in Fig. 2. 7. Using 
the symmetry of the unit cell a 2D model can be presented. 



Fig 2.7.Unit cell and axisymmetrical FE model 


A penny-shaped crack resulting from fiber-end debonding is assumed to exist. The mechanical loading is 
tension. The matrix is described as an elastic -plastic material. Under loading micro-void nucleation, growth and 
coalescence are observed. Near the corner of fiber a local plastic zone is formed. If the plastic strain reaches the 
critical fracture strain value crack will propagate. The value of the critical strain can be determined 
approximately. The angle of crack affects the direction of crack propagation; the angle is described by the 
location of the maximum stress. The paper presents mechanical loading, thermal loading and the combination of 
the two above. Results show those thermal residual stresses and strains strongly affect the condition of the 
interface and change the location of the maximum equivalent strain. 

Laminated structures also require investigating the reason of matrix damage. A [0°, 90°3, 0°] T cross-ply 
structure is investigated in [19], in the case of epoxy matrix and glass fiber. The cross-ply contains long 
continuous fibers. The FE model is shown in Fig. 2. 8. It is assumed that fibers are uniformly distributed in the 
matrix and have the same radii. The loading is strain controlled; the boundary conditions are also depicted in 
Fig. 2. 8. Attention is focused on the interior domain; no free edge effect is included. 

Fibers are linear elastic; matrix behaves as nonlinear viscoelastic material. Both constituents are homogeneous 
and isotropic. 


A 


TUT K>.TGI>: 


• rFDT.'Trr 


K , r.'.'.YvVT' 



ri* * V * t jN 

\ i«-»\ 0 ‘-U 

* •; o 


• V V ■ i > V Tv \ . .. X \ ! J : .* 

\ v r k / // 


s s 


-'W 


Fig 2.8 3D FE model and location of the damage, [0°, 90°3,0°]t 

Criterion for the epoxy matrix describes if the maximum principal strain reaches the critical strain value, damage 
progression will occur in matrix. FE results indicate that local damage zones arise in two areas, shown in Fig. 2. 8. 
Increasing loading both damage areas grow and coalesce. Applying linear analysis the damage occurs at higher 
loading. So that the nonlinear model describe more accurately the real process. Experimental datas support this 
fact. Nonlinearity plays important role in the analysis of composites. 
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2.3.4 Modelling of composite structures with holes and notches 

Numerous applications of composites require the presence of holes and cutouts. These can cause inhomogeneous 
stress distribution in the structure. Cracks and flaws act important role, because these affect the strength of the 
composite. Crack can be caused by the manufacture process or loading. Concerning failure criteria, three basic 
types of approaches can be found in literature [8, 10]. 

Fracture mechanics models: 

These models assume a localized damage emanating from the hole, represented by an equivalent crack, which 
causes failure when it reaches a critical size. The strength of the laminate is related to the fracture toughness or 
the strain energy release rate. The fracture toughness and the critical size seemed to be hole-size and material 
dependent. This is the basic idea of LEFM (Linear Elastic Fracture Mechanics). 

Point and average stress criterion: 

It is assumed that failure occurs when the stress in one point at a certain characteristic distance or the average 
stress along the characteristic distance from the hole attains the unnotched strength of the laminate. Both 
criterions can be expressed in closed-form equations. Experimental data shown those characteristic distances are 
also hole-size and material dependent. 

Progressive damage models: 

Throughout the loading process damage and resulting changes in the stress distribution occur. The appropriate 
combinations of failure modes, failure criterion and property degradation law are required. 

Open-hole tensile strength of quasi-isotropic laminates is investigated by MORAIS in [8]. Stress analysis near 
edges and discontinuities can only be accurately performed using 3D models. 2D models give non -negligible 
errors in laminates. The central hole has the effect of stress concentration, the stress concentrates around the hole. 
Considering this fact MORAIS presented a mixed 2D/3D model. The model is shown in Fig. 2. 9. The specimen 
was modeled with 2D elements, except for a narrow zone around the hole, where 3D brick elements were used. 
The connection between the 3D and 2D region was done imposing identical in-plane displacements on the 
boundary nodes of identical in-plane coordinates and using interface elements. Longitudinal strength of 
laminated structures depends on the size. Failure criterion does not follow this fact. First a [0745 7907-45 °]s 
laminate was investigated. The loading and boundary conditions are depicted in Fig. 2. 9. The laminate is 
transversely isotropic; properties can be obtained using classical laminate theory. The inhomogeneous stress 
distribution can be 



obtained along a central ligament. The 3D stresses considered are layer average stresses. The stress concentration 
depends on the hole-size. Comparing the 2D and the mixed models the mixed model shows higher stress 
concentration values near the hole than the 2D model. 2D models under predict stress concentrations, further 
from the hole approximately identical values were obtained by both models. The strength of structure with 
central hole depends on the longitudinal strength of the laminate. In situ fiber strength can be estimated using 
experimental datas and the mixed model. 

Laminates with triple parallel-arranged and equally spaced notches can be found in [20]. Glass- and carbon- 
fiber/epoxy composites containing through or surface notch were investigated. Fig. 2. 10 shows the laminate and 
the FE model. The FE model is only quarter of the structure, the stacking sequence is [07907079070°]s. The 
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laminated composite behaves as orthotropic linear elastic continuum, the interlaminar matrix is isotropic. The 
paper presents results taking account of the interlaminar matrix and not respectively. Notch can be surface or 
through notch. The mesh was gradually refined to the notches. Stress concentration factors (SCF) were 
introduced for the central and the deviate notch varying the distance between the two notches. The depth of layer 
describes the surface notch. Results show that stress concentration factors near the notch tips change from high to 
low then revert to high with spacing between both neighboring notches. 



\:TJ*w4L_> 


Fig 2.10 Laminated structures and FE model with triple parallel notches 

Stress concentration factor always larger at the deviate notch. Depth of the notch increases the SCF. In order to 
validate the FE predictions a total of twenty specimens with triple parallel notches were examined. Failure tests 
show the change of the failure load exhibits an opposite trend with the notch stress concentration factors. Tests 
also show that failure always occurs at one of the two deviate notches. 

A notched laminated structure, including splitting and delamination was investigated by WISNOM and CHANG 
[18]. A cross ply laminate of carbon/epoxy with lay-up [9070°] with a center crack shown that splitting occurs in 
the 0° plies perpendicular to the notch together with narrow areas of delamination. The [9070°] laminate was 
modelled with four nodded plane elements. The elements are plane stress ones. Duplicate nodes were used on 
either side of any ply interfaces where delamination is expected. The same interface element was also used to 
model splitting in the 0° ply. The schematic illustration and the FE model are depicted in Fig.2.11. Separate 
elements were used to represent the 0° and 90° plies. Duplicate coincident nodes were used for the two plies to 
allow relative displacement to occur during delamination. These nodes were connected with interface elements. 
Initially linear elasticity was assumed, but it became apparent that the nonlinear in -plane shear response was very 
important. The interface was modelled with nonlinear springs between the coincident nodes. The springs are 
initially assumed to be elastic. When a critical relative displacement between plies is reached the interface is 
assumed to fail. 

The tensile loading was applied in increments. Initial results show splitting initiated from the end of the notch, 
and then progressively increased in length. 



FE MODEL 

Fig.2.11 Schematic illustration of damage and FE model 
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It was accompanied by the formation of a narrow zone of delamination. The extent of delamination is slightly 
less than the split length. It is the splitting that drives the delamination. The analysis was repeated with refined 
FE mesh. Results show that the analysis is not sensitive to the mesh size. The difference between linear and 
nonlinear analisys was also investigated. Using linear analysis gives unreal results; hence non-linearity acts 
important role. The analysis was also repeated with the spring stiffness halving to examine the effects of 
interface parameters. In this case the damage progression is similar to the original case, but the size of the yield 
zone has greatly increased. 

2.3.5 Interlaminar delamination and transverse cracking 

Damage models can be constructed using proper combination of the global element types. NISHIWAKI et al. 
developed quasi-three-dimensional model for the modelling of interlaminar delamination and transverse cracking 
independently [9]. The FE model is depicted in Fig. 2. 12 and consists of shell and beam elements, which 
represents fiber and resin respectively. The fiber reinforcement is concentrated in the shell elements with the 
maximum fiber- volume fraction of 90, 7%. Beam elements represent the interlaminar matrix. The loading is 
lateral compression. 



Fig 2.12 Four layered quasi-three-dimensional cylinder model 

The cylinder model was investigated with different stacking sequences. Results show different loading limits, 
transverse cracking occurs if a given shell element reaches the yield stress and delamination can be observed if one of 
the beam elements yields. Increasing compressive load more and more shell and beam element reaches the yield stress. 
Interlaminar delamination is caused by shear stresses. Finally the estimation of the critical strength was performed. The 
critical strength was the load level corresponding to initial yielding of beam elements along the axial direction 

Interlaminar fracture is one of the major problems for fiber composites. Its occurrence greatly reduces the 
stiffness of a structure. Modelling of interlaminar crack can be performed 
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Fig2.13*End notched flexure specimen and FE mesh 


using end-notched flexure specimen (ENF) under in -plane shear deformation mode [17]. Applying the ENF the 
fracture toughness can be measured. The aim is to provide quantitative information on interlaminar crack growth 
in fiber composites and understand the role of fibers on the stress distribution in the surrounding matrix. The 
used FE model and the meshed structure are shown in Fig.2.13. The model contains only a single fiber on each 
side of the starting defect. The main concern is the role of the blunt starting defect, introduced by the notch. The 
in-plane shear loading was introduced by assigning a fixed displacement to the top and bottom surfaces. 
Boundary conditions are also shown in Fig.2.13. Experimental ENF tests were conducted to verify FE results 
using two different specimen. The maximum principal stress shows that the highest stress occurs around the 
corners of the starting defect, extending towards the interface. The stress distribution indicates two mechanisms 
that are responsible for the fracture initiation. Fig.2.14 depicts that the fracture can be initiated either from the 
corner of the starting defect or at the interface. The latter occurs when the interfacial strength is lower than a 
critical percentage of the matrix strength. This percentage can be obtained from the stress distribution around the 
notch. 
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Fig.2.14 Local damage process 

When the crack is initiated from the interface, its growth is expected in both directions. The further crack can 
grow towards the starting defect and interfacial debonding occurs. If the interface strength is sufficiently high, 
the crack initiation is caused by tensile failure of the matrix around the starting defect. The crack propagates 
towards the interface because of the shear failure. It is expected that the interface strength can significantly 
influence the measured strain energy release rate. 


2.3.6 Modelling of fiber-fragmentation 

A finite element analysis was carried out by DAVIS and QU to simulate the fiber fragmentation process in a 
single-fiber composite (SFC) specimen [3]. A single fiber is embedded in a uniform matrix material. Under 
monotonic axial loading the fiber will fracture. As the load increases further, more fiber fracture takes place as a 
result of load transfer from the matrix to the fiber. The fiber-fragmentation process will stop because either 
interfacial debonding or the fragments are so short that not enough loads can be transferred to the fiber for 
further fracture. This saturated length is related to the interfacial strength. The SFC fragmentation test is suitable 
to characterize the in situ interfacial strength and the understanding of the micromechanics in SFCs can be 
transferred to real composites with higher fiber-volume fractions. The A1 matrix is assumed to be elastic/plastic; 
the interface is modeled by elastic/perfectly plastic spring layer. Fig.2.15 illustrates the SFC and the applied 
undamaged and damaged FE models. As the applied load is increased, the fiber will fracture when the fracture 
strain is reached. The fiber fracture creates a penny-shaped crack, creating a localized plastic zone in the matrix. 
The number of fracture per unit length will increase; i.e. the fiber-fragment length will get shorter and shorter as 
load increases. The shear stress on fragments, which depends on the interface and the fiber strength. The SFC 
specimen containing many fractures can be simplified into a unit cell model containing one fragment. Only one 
quarter of the unit cell needs to be modelled for finite element analysis. A series of finite element models of the 
unit cells with decreasing length was used to analyze the fragmentation process. To simulate the successive fiber 
fracture, a series of 21 SFC unit cells with decreasing length were used, applying solid elements. The procedure 
was the following. 

The fiber strength was calculated and the load was increased until the axial stress in the middle of the fiber 
fragment exceeded the fiber strength. With the next specimen the calculation was repeated. This was done 
without and with fiber fracture on the model. It is found that although the radial and tangential stresses in the 
matrix vary radially, the magnitude of the longitudinal stress is much higher. Therefore the latter dominates the 
von Mises stress. One notable difference between the unfractured and fractured models is that the longitudinal 
shear stress was not zero in the latter case. 



Fig.2.15.Finite element modeling of fiber fracture 

It seems that the stress-strain curve is not a good indicator of the interfacial shear strength. The applied effective 
strain/normalized fiber-fragment length is much suitable to describe the effect of interface strength. Finally it is 
possible to estimate the interface strength by inspecting the fiber-fragment lengths from a SFC test at a given 
effective strain. 
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III. FINITE ELEMENT MODELLING 


k 


t 


Fig 3.1 Representative volume element method 

3-Dimentional Finite Element models are developed with governing boundary conditions to study the response 
of the unit cell due to the external loads. The Finite Element Software ANSYS 10.0 is successfully executed for 
the analysis. The present finite element model is validated by comparing with the exact mechanics of materials 
results for the case of perfect bonding; the analysis is extended to predict the mechanical properties of the 
composite for different volume fractions with perfectly bonding at fiber-matrix interface. 

The 1-2-3 coordinate system shown in Fig. 3. 2 is used to study the behaviour of a unit cell (The direction 1 is 
along the fiber axis and normal to the plane of the 2D plane given in Figs.3.2&3.3). The isolated unit cell 
behaves as a part of a larger array of unit cells. 

It is assumed that the geometry, material and loading of the unit cell are symmetrical with respect to 1-2-3 
coordinate system. Therefore, a one fourth portion of the unit cell is modeled (Fig. 3.3) for the prediction of 
mechanical properties. The 3D Finite Element mesh on one fourth portion of the unit cell is shown in Fig. 3.3. 




t S y 




Figure 3.3 3-D Finite element mesh 


3.1.3. Element Type 

The element SOLID95 of ANSYS 10.0 used for the present analysis is based on a general 3D state of stress and 
is suited for modeling 3D solid structure under 3D loading. The element has 20 nodes with three degrees of 
freedom per node (UX, UY and UZ). 

3.1.4 Materials 

The T300 fiber and polymer matrix materials with following properties are used. 

• T300 Fiber: Ei = 220.632 GPa, E 2 (= Es) = 13.789 GPa, V12 (= vis) = 0.2, V23 =0.25, G 12 (= Gis) = 8.963 GPa, 
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G 23 = 4.826 Gpa, ai = -0.99x10-6 / K, 
a 2 (= as) = 10.08x10-6 /K 

• Polymer Matrix: E = 5.171 GPa, v = 0.35, a = 72x10-6 /K 

Boundary Conditions: Due to the symmetry of the problem, the following symmetric boundary conditions are 
used. 

At x = 0, Ux = 0; at y = 0, U y = 0; at z = 0, Uz = 0 

In addition, multipoint constraints are given so that the plane faces of the unit cell remain plane after 
deformation. 

The Young’s modulii are determined using the equations 

El = Ol / 8i; E2 = G2 / 82; E3 = G3 / S3 

The Poisson’s ratios are determined using the equation 

V12 = - 82 / Si; V13 = - S3 / Si; V21 = - SI / S2; V23= -S3 / 82, V31 = - Si / S3; V32 = -82 /S3 

3.2.3. Finite Element Model 

The 1-2-3 Coordinate system shown in Fig. 3. 5 is used to study the behavior of unit cell where 1- is the fiber 
direction, perpendicular to the plane of the figure. The isolated unit cell behaves as a part of large array of unit 
cells by satisfying the conditions that the boundaries of the isolated unit cell remain plane. 

3 * 



Fig 3.4. Isolated Unit Cell 


It is assumed that the geometry, material and loading of unit cell are symmetric with respect to 1 -2-3 coordinate 
system. Therefore, a one-eighth portion (half in cross section and half along the length direction) of the unit cell 
is modeled for the analysis (3.6.). 

* >'<-> 



3.2.4. Geometry 

The dimensions of the finite element model are taken as x = 100 units (2-in-plane transverse direction), y = 200 
units (3-out-of-plane transverse direction) and z = 10 units (1 -Fiber direction). The radius of fiber is varied 
corresponding to the volume fraction. For example, the radius of the fiber is calculated as 6F8 units, so that the 
fiber volume fraction becomes 0.30 

3.2.5. Element Type 

The element used for the present analysis is SOFID 95 of ANSYS 10.0 software, which is developed, based on 
three-dimensional elasticity theory and is defined by 20 nodes having three degrees of freedom at each node: 
translation in the x, y and z directions. 
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3.2.6. Materials 

The top fiber in the unit cell is taken as boron and the bottom fiber as graphite. The properties of the constituent 
materials are given in Table3.1, which is shown here under 

Table 3.1: Properties of Constituents 
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3.2.7. Loading 

Uniform temperature load of IK is applied on the entire volume of the FE model. 

3.2.8. Boundary conditions 

Due to the symmetry of the problem, the following symmetric boundary conditions are used. 

At x = 0, Ux = 0; At y = 0, U y = 0; At z = 0, Uz = 0 

In addition, the multipoint constraints are imposed on the boundaries of the unit cell to make them plane after 
loading. 

IV. RESULTS AND DISCUSSIONS 

The predictions of elastic constants are carried out by FEA and mechanics of materials approach. The results of 
theory of mechanics of materials and FEA are tabulated in Table 4.1, 4.2 for perfectly bonding case. It is 
observed that 3 decimal accuracy i.e. found to a close agreement. 


Table 4.1: Comparison of El for perfect bonding 



Analytical 

solution 

+ 
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technique (G 
Pa) 
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error 
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Table 4.2: Comparison of E2 for perfect bonding 
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(G Pa) 



0.15 

86.125 

87.214 

1.24 

0.3 

110.240 

111.098 

0.77 

0.45 

144.001 

143.669 

0.02 

0.6 

187.408 

187.318 

0.04 

0.75 

244.000 

244.834 

0.34 


4.1FINITE ELEMENT RESULTS FOR T300-EPOXY AT VOLUME FRACTIONS 0.6, 0.7, 0.8, and 0.9: 

The 2D and 3D finite element models and the deformed and undeformed shape of finite element model and 
Stress and strain distribution at various volume fractions are shown from figures 4.1.1 to 4.1.27. 

For volume fraction 0.6 


AN 



Fig.4.1 3D finite element model 



Fig 4.3 Strain distribution 




Fig: -4 .5 3 O finite element .model. 
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Fig: 4.8 Stress distrfbutijQ^ 



Fig: 4.11 Strain on YZ plane 


For volume fraction 0.8 



Fig: 4.1.15 2Dfinite element model 



Fig: 4.9 Strain distribution 

VWWSAAAAAAAA#VWWV 




Fig: 4.13 Deformed and un deformed shape 
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Fig: 4.16 Stress distribution 
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n g: 4.19 3D finite element model 
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Fig: 4.20 Deformed and undefoi^ed shape 






Fig: 4.21 Stress dUXlJJfeJitlaBt 


For volume fraction 0.9: 



Fig: 4. IS 2D finite element model 



www.theijes.com 


The IJES 


Page 69 


Micro Mechanical Modeling Of Fiber / Epoxy... 



V.FINITE ELEMENT RESULTS FOR HYBRID COMPOSITE LAMINA SUBJECTED TO 

THERMAL LOADING: 

Sufficient number of convergence tests is made and the present finite element model is validated by comparing 
the coefficients of thermal expansion of graphite -epoxy lamina predicted using FEM with the value obtained 
from Alternate Rule of Mixtures and found close agreement (Figs. 5.1 and 5.2) Equations for ai and ai ARM is 
given as 

oil =[[ai f Ei f -a m E m ] V f+ a m E m ] / [[ Ei 1 - I 

m ]Vf+ E m ] 

m + , f m \ -v r „ //T n f, ,m m ^ . f. , _ 

a2 = a (ct2 - a ) Vf + ((Ei V - E Vi ) / Ei 

(a m - a i f ) ( 1 - Vf) Vf 
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Fig 5.2 Validation for Pi with Alternate Rule of Mixtures 
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Fig. 5.1 shows the variation of the coefficient of thermal expansion of the composite in the fiber direction (ai). It 
is observed that ai decreases with Vf for all the three composites. This is because the coefficient of thermal 
expansion of matrix material is more than the fibers. It is also observed that the values of ai are more for the 
boron epoxy followed by hybrid-epoxy than graphite -epoxy at all volume fractions. The reason for this is the 
values of ai are positive for boron and epoxy, and negative for graphite epoxy. The variations of transverse 
coefficients of thermal expansion with respect to Vf are shown in Figs. 5.2 and 5.4. In these cases also the value 
of coefficient of thermal expansion decreases with increase in volume fraction. The values of ai and a3 are 
observed to be more for graphite -epoxy followed by hybrid-epoxy than boron-epoxy as the transverse coefficient 
of thermal expansion for graphite material is more compared to boron material. 



Volume Fraction (Vf) 


—♦—BORON EPOXY 
-■-GRAPHITE EPOXY 
-A- HYBRID EPOXY 


Fig 5.3 Variation of Cl: 'with respect to volume fraction 


Fig 5.5 Variation of □» with respect to volume fraction 


Validation of thermal expansion at volume fractions 0.15, 0.3, 0.45: 

The 2D and 3D finite element models and the deformed and undeformed shape of finite element model and 
Stress and strain distribution at various volume fractions are shown from figures 5.6 to 5.20. 

For volume fraction 0.15 
Fig: 5.6 2D finite element model 
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Fig: 5.7 3D finite element model 
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Fig: 5.8 Deformed and undeformed shape 



Fig: 5.9 Stress distribution 



THERMAL STRESSES IN A SINGLE ELASTIC FIBER EMBEDDED IN AN 
INFINITE MATRIX 

The comparison is based on the normal and shearing stress distribution along the fiber. The analytical results 
come directly from Eqs (3.14) and (3.15). The numerical results were obtained using the finite element analysis 
results from ANSYS. Different cases were considered as a result of varying the material properties, ratio b = Em/ 
Ef, and the radial dimensions of the fiber, ratio 1 = r 0 /l. Figures 4.3.1. to 4.3.6 compared the numerical and 
analytical results of the normalized residual axial and interfacial shear stresses. 

It is observed that the numerical results showed the same qualitative behavior of the analytical solution, serving 
as a validation tool against lack of experimental results. There was a large discrepancy in the results as the ratio z 
/I approached one. This is attributed to two reasons. The first one is that the analytical solution considered the 
axial stress at z /I = 1 small enough to be assumed zero; however, in the numerical analysis is observed that the 
axial stress at z /I = 1 is very small but not necessary zero. Therefore, there exists a small discrepancy in the 
stresses between the numerical and analytical solution mainly when z/ 1 approaches one. This fact is especially 
seen in the case of short fiber aspect ratio. The second reason for the discrepancy in the results at the end is due 
to a singularity. Singularities occur at locations where a drastic transition in geometry or a difference in the 
properties of the materials is present. The best way to deal with this is by increasing the number of elements 
particularly in the affected region or by using singular elements in the finite element model. The latter approach 
has been used in fracture mechanics to determine the stress intensity factor at the crack tip. The approach used in 
this work was to increase the mesh size in the fiber and in the matrix, especially at the area close to the interface, 
so the discretization was maximized not only at the tip but along the whole interface. Then the results at the tip 
were simply discarded, knowing that the erratic behavior was already expected. This erratic behavior was 
observed around the last three percent of the fiber’s length, close to z/ 1 = 1, and was particularly noticed in the 
shear stress numerical results rather than in the axial stress. This was also expected since the shear stress was 
obtained from a single nodal value at the interface, while the axial stress was the result of the average across the 
fiber radius. 


15 
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region is not actually known. For the case when multiple fibers are present, the ratio of the perturbed radius to 

fraction. Still the optimum value for the distance between 
fi bers is unknown. Equation s (3. 17) and ( 3.18) have p roven satisfactory in approximating the perturbed region in 
this work and in similar type of problems. They provided a good qualitative solution to the residual embedded in 
a soft material. This stresses when compared to the numerical solution from ANSYS. This applied to all 
conditions; short or long fiber embedded in a soft matrix, and short or long fiber embedded in a stiff matrix. 
However, the best agreement between the analytical and numerical solutions occurred for the case modeling a 
long fiber validated the original assumption of taking X« 1 and p« 1. The error for the problem when 1 = 0.005 
and b = 0.001 started at around two percent in the axial stress at z /I = 0, stayed within five percent through z /I = 
0.6, then increased to ten percent at z /I = 0.75 , twenty percent at z/1 = 0.85 and continued increasing to the end 
of the fiber. The shear stress started around twenty percent at z /I = 0, continued decreasing to less than one 
percent at z/1 = 0.9, then increased again to the end. The reason for the discrepancies in the stresses as the end of 
the fiber is reached was discussed earlier. 



Fig 5.21 Comparison of normalized axial stress for 0 = 0.005, 0=0.01, B Dm =-100e-6 1/deg centigrade 

□ T=600 deg cen 

Fig 5.22 Comparison of normalized shear stress for 0=0.005, 0=0.01, 0 Dm =-100e-6 1/deg centigrade 
□ T=600deg cen 



Fig 5.23 Comparison of normalized shear stress for 0=0.005, 0=0.001, 0 Dm =-100e-6 1/deg centigrade 
0 T=600 deg cen 
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Fig 5.24 Comparison of normalized shear stress for □=0.005, □=0.001, 0 HU =-100e-6 1/deg centigrade 
0 T=600 degcen 



Fig 5.25 Comparison of normalized shear stress for 0=0.01, |3=0.001, 

□ T=600 deg cen 


B Dm =-100e-6 1/deg centigrade 



Fig 5.26 Comparison of normalized shear stress for 
□ =0.01, 0=0.001, 0 Dm =-100e-6 1/deg centigrade OT=600degcen 
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Ansys results for Stress for □=0.005, D=0.01, 0 Dm = - 100e-6/deg c FOR Z=0.2, 0.3, and 0.9: 

The 2D finite element model and the deformed and undeformed shape of finite element model and Stress and 
strain distribution at various Z values are shown from figures 5.7 to 5.18. 


For Z=0.1 


CONCLUSIONS AND FUTURE SCOPE 


The following conclusions are drawn from the present work 

1 . Engineering properties of T300 - epoxy fiber reinforced composite lamina have been predicted 
using mechanics of materials approach based finite element method applied to representative volume elements in 
the form of square unit cells with fiber matrix perfect bonding. The average percentage error in between 
Analytical and Finite element method is given as 0.01. 

2. Micromechanical analysis of a hybrid fiber-reinforced composite due to thermal loading is studied 
in the present work. The coefficients of thermal expansion of the boron-graphite-epoxy hybrid composite are 
predicted using three-dimensional finite element method. The coefficients of thermal expansion are found to be 
decreasing with increase of volume fraction (Vf) as the value of coefficient of thermal expansion is more for the 
matrix material. 

3. The numerical analysis confirmed the approximate analytical solution to the thermal stress along a 
fiber embedded in a matrix as presented in this work. It is observed that the accuracy of coincidence of both 
results depends on the dimensionless ratios X=r a /l and P = Em/Ef. 

However, both approaches coincide qualitatively and quantitatively, even though substantial discrepancies might 
appear at particular points mainly due to the strong singularity present on this type of problem. The percentage of 
error decreased as p decreased typically for longer fibers. This was 

especially observed for the interfacial shear stress error. As P decreased from 0.01 to 0.001 the maximum error 
decreased from fifty to the low twenty percent for the case when X = 0.005. This agreed with the initial 
assumption of using e much less than unity are for the case when X = 0.01 and p 
= 0.001. 

Future scope: 

1 . The influences of fiber-matrix interface debond on 

the behaviour of fiber reinforced composite lamina and considering the matrix material as elastoplastic, 
viscoplastic, elastoviscoplastic is the topic for future work. 

2. The composite is built ply’s with different orientations are considered as future scope. 

3. The next step would be to study the interaction between two or more fibers, finding an expression for 
the optimal distance of the perturbed region that produces the best representation of the stress transfer between 
them. These topics will be discussed in future work. 
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